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FOREWORD 


This is a progress report on the research project “Surface Modeling and Optimization Studies of 
Aerodynamic Configurations”. Within the guidelines of the project, special attention was directed toward 
research activities in the area of “Aerodynamic Shape Optimization of a HSCT Type Configuration With 
Improved Surface Definition”. The period of performance of this specific research was May 1, 1993 
through June 30, 1994. 

This work was supported by the NASA Langley Research Center through Cooperate Agreement 
NCC1-68. The cooperate agreement was monitored by Dr. Robert E. Smith Jr. of Analysis and 
Computation Division, NASA Langley Research Center, Mail Stop 125. 
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ABSTRACT 


AERODYNAMIC SHAPE OPTIMIZATION OF A HSCT TYPE 
CONFIGURATION WITH IMPROVED SURFACE DEFINITION 

Almuttil M. Thomas 1 
Surendra N. Tiwari 5 

Two distinct parametrization procedures of generating free-form surfaces to represent aerospace- 
vehicles are presented. The first procedure is the representation using spline functions such as Non- 
Uniform Rational B-Splines (NURBS) and the second is a novel (geometrical) parametrization using 
solutions to a suitably chosen partial differential equation. The main idea is to develop a surface which 
is more versatile and can be used in an optimization process. Unstructured volume grid is generated 
by an advancing front algorithm and solutions obtained using an Euler solver. Grid sensitivity with 
respect to surface design parameters and aerodynamic sensitivity coefficients based on potential flow is 
obtained using an Automatic Differentiator precompiler software tool. Aerodynamic shape optimization 
of a complete aircraft with twenty four design variables is performed. High speed civil transport aircraft 
(HSCT) configurations are targeted to demonstrate the process. 


t Graduate Research Assistant, Department of Mechanical Engineering, Old Dominion University, Norfolk, VA 23529-0247. 
§ Eminent Professor, Department of Mechanical Engineering, Old Dominion University, Norfolk, VA 23529-0247. 


iii 



TABLE OF CONTENTS 


E2S£ 

FOREWORD il 

ABSTRACT iii 

1. INTRODUCTION 1 

1.1 Engineering Design and Surface Parametrization 2 

2. SURFACE MODELLING AND GRID GENERATION 4 

2.1 Rational B-Spline Curves and Surfaces 4 

2.2 The PDE Method 6 

2.3 Graphic Interface 9 

2.4 Unstructured Grid Generation and Method of Solution 9 

3. GOVERNING EQUATIONS FOR POTENTIAL FLOW SOLUTION 17 

4. COMPUTATION OF SENSITIVITY TERMS WITH RESPECT TO DESIGN PARAMETERS . 19 

4.1 Grid Sensitivity with respect to Design Parameters 20 

4.2 Optimization Problem 21 

5. CONCLUSIONS 22 

6. REFERENCES 26 


IV 



1. INTRODUCTION 


The design of aerospace vehicles requires the solution of mathematical models to represent surfaces 
which are systems of either algebraic, differential or integral equations [1]. The defining characteristic of 
these surfaces is that they must be changeable by the designer. The use of rational polynomial functions 
for representing curves and surfaces in CAD/CAM applications is becoming increasingly important 

The type of polynomial functions used gives various forms of curves and surface representations, for 
example, coons patch, Bezier curves and surfaces, B-splines, and rational B-splines. Free-form surface 
design using such curve and surface representations has been discussed by a number of authors. The 
method of using rational B-spline curve to produce a variety of different surface shapes is described by 
Tiller [3-4]. Woodard [5] presents a variety of techniques for producing free form surfaces and goes 
on to describe, in some detail, interactive skinning technique using B-spline surface interpolation. The 
increasing popularity of rational Bezier and B-spline forms is due to the fact that they offer one common 
mathematical form for precise representation of standard analytical shapes such as lines, conics, circles, 
planes and quadric surfaces [6]. 

Free form surfaces can also be generated as a solution to a suitably chosen Partial Differential 
Equation (PDE), and this has been investigated by Bloor and Wilson [7-9]. The method uses PDEs as a 
means of producing blend surfaces. By regarding blend as a solution to a boundary-value problem and 
by choosing appropriate boundary conditions, it was demonstrated that a solution to an elliptic PDE gave 
a smooth blended surface that has the required degree of continuity with the surface to which it joined. 
The surfaces generated by this method are expressed parametrically, often in terms of transcendental 
functions of the surface parameters rather than simple polynomial expressions, and hence the resulting 
surfaces tend to be smooth. 
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Both of the above described methods have been successfully applied to various types of surface 
shapes but there has been relatively little work published on the use of these forms for design of aerospace 
vehicles. The main goal of this study is to focus on developing surfaces for aerospace vehicles using 
NURBS and PDE which can be later used to compute computational grids in a very rapid manner. 
The ultimate goal is to produce surface and grids suitable to optimize designs in conjunction with the 
numerical simulation of the physics. 

Grids are generated to discretize the solution domains of the physical-mathematical models so that 
numerical solutions can be obtained. A grid is defined as a set of points with appropriate connections 
between the points. The points act as reference positions within the field at which the flow variables are to 
be computed and the connections between the points act as pathways for transferring information around 
the computational domain. In a structured grid the connectivity between the points is implicitly defined 
through a curvilinear coordinate system. In an unstructured grid, the connectivity is arbitrary and therefore 
must be explicitly specified. Solution methods that utilize a structured grid are generally more efficient 
than methods that utilize an unstructured grid. However unstructured grids provide a much greater degree 
of flexibility than is available with a structured grid. In particular, unstructured grids can discretize a 
highly complex domain easily and are suitable for performing localized grid enrichment for solution 
adaptation. In this study unstructured grids around various surfaces considered are generated using the 
advancing front technique [10]. This method was selected because it does not require a separate library 
of modules to distribute grid points throughout the domain in advance like the Delauny triangulation. 

Discussion on engineering design, surface modelling and grid generation technique are briefly 
presented. This is followed by governing equations for potential flow solution and finally calculation of 
sensitivity terms and optimization problem is presented and discussed. 

1.1 Engineering Design and Surface Parametrization 

In developing mathematical models and numerical solution techniques, parameters that characterize 
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the discipline have evolved [11]. The parameters are divided into independent parameters Pi„d and 
dependent parameters P dep . Independent parameters characterize discipline physics and the solution 
domain geometry. Examples of independent parameter for aerospace vehicle are: Mach number, Reynolds 
number, wing sweep, mean camber, maximum wing thickness, fuel weight, chord length, and panel 
thickness. Independent parameters broadcast information to specify conditions in the solution domain. In 
the case of geometry, independent parameters define the grid and vehicle surface. Since a mathematical 
model in one discipline can require input from other disciplines, this input may be classified under physics 
parameter. Dependent parameters are usually integration of dependent variables in mathematical models. 
Examples are lift, drag, weight, and wing volume. 

For an aerospace vehicle such as High-Speed Civil Transport (HSCT) the traditional approach to 
design is for aerodynamics and performance disciplines to initially create the vehicle surface [12]. The 
process is to define the planform, wing, fuselage, engine nacelles, and major control surfaces with 
aero/performance independent-design parameters. Approximately 50-100 independent parameters are 
required to specify a rough vehicle surface design. Usually a sparse set of points on component surfaces 
which can be thought of as a coarse grid becomes the surface description for analyses. In the CAD 
process a considerable amount of information is added which is not provided by the low-level analyses. 
These additions are in the form of point movements and point creations to achieve desirable surface 
characteristics such as smoothness and blending of one component into another. These two aspects are 
kept in mind in developing surfaces for the different geometries considered here. 
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2. SURFACE MODELLING AND GRID GENERATION 


Two different approaches have evolved in the development of parametric curves and surfaces. They 
are referred to here as “interpolative” and “approximative”. In an interpolative representation, points 
and derivatives on the curve or surface are used to control the formula defining the curve and surface. 
Lagrangian and Hermite interpolation formulas are examples of this approach. In an approximative 
approach, points not necessarily on the curve or surface control the formula defining the curve or surface. 
Bezier and B-spline representations are examples of this approach. In the design process using an 
interactive CAD system, the approximative approach is highly advantageous. The basic formulation of 
this approach is presented here. 


2.1 Rational B-Spline Curves and Surfaces 


A NURBS curve [4] is a vector valued piecewise rational polynomial fimction of the form 

£ oj t V,B t ' k (u) 

Q ( «) = ^ ( 2 - 1 ) 

. £ UiB itk (u) 

i=0 

where w, are called the weights, the V, are the control points and B i k (u) are the normalized B-spline 
basis functions of degree k defined as 




if UjKuKU'+j 
otherwise 


( 2 . 2 ) 


5, •,*(«) = — — —Bi, k-i(u) 


+ 


u i+k u i 
u i+k+l “ u 


“ Ui+i 


B i+ i, P -i(u) 


where Ui are the so-called non-uniform knots forming a knot vector. 


(2.3) 
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An interactive program based on the formulation given by Eqs. (2.1) - (2.3) has been developed. 
In this program, after prescribing an initial set of control points, the designer can pick and drag these 
points and simultaneously observe the change in the shape of the curve. Weights can also be specified at 
each control points and the curve can be modified. As a simple example, this approach is applied to the 
design of cambered airfoils. The airfoil is controlled directly with NURBS control points and weights. 
Figure 2.1 shows a six control point definition of the cambered airfoil. The points at the leading edge and 
trailing edge are fixed. Two control points at 0% chord are used to affect the bluntness of the section. 
The effect of the movement of the control points to create another airfoil is shown in Fig. 2.2. Figure 2.3 
shows the effect of increasing the weight of the middle control point. It is seen that the curve is pulled 
towards the control point. An arc length based discretization of the unit line is used for the knot vector. 

A NURBS surface [6] is the rational generalization of the tensor product nonrational B-spline surface 
and is defined as 




Of \ *=0]-0 

S(u,v) = — - 


(2.4) 


\ — 0 j = 0 


where are the weights, Vy form a control net, and Bi,k(u) and Bj <q (v) are the normalized B-splines 
of degree k and q in the u and v directions, respectively. The knot vectors are 

U = {0,0,...0,u k +i,...,u r -k-i,l, 1,-1} 


V = {0,0,...0,tt, +i , 

where the end knots are repeated with multiplicities k+1 and q+1 and r = n + k + 1 and s = m + q + 1. 

n m 

A NURBS surface has the property Y E v) = 1 and reverts to a B-spline when all 

i=0 j=0 

the weights are 1. It has the advantage of being able to represent free form surfaces, and with the proper 
choice of weights, conic surfaces. Surface skinning technique [5] is used to obtain the NURBS surface. 
The task of skinning is to fit a surface through an ordered set of space curves, called as section curves. 
The positioning of section curves in the three-dimensional space is customarily done with respect to a 


5 



spline curve, from which appropriate orientation vectors can be automatically computed. The surface 
skinning technique is used to define NURBS surface. An ONERA M6 wing is used for this case. The 
wing has a leading edge sweep of 30 degrees, an aspect ratio of 3.8, taper ratio of 0.56, and symmetrical 
airfoil sections. The wing defined by three wing sections and nine control points per section is shown in 
Fig. 2.4. The NURBS surface generated using these control points is shown in Fig. 2.5. 


2.2 The PDE Method 


The PDE method generates a surface X_ in Euclidean 3-space, which is a function of two parameters, . 
i.e., X_ = (x(u,v), y(u,v), z(u,v)). The surface is obtained by solving a partial differential equation (PDE), 
in parameter u,v space, subject to boundary condition on X_ and its normal derivative with respect to 
u and v. In general the order of PDE determines the number of derivatives of the unknown function 
that must be specified in the boundary condition. If control over both shapes of the curves bounding 
the PDE surface patch and the directions and magnitude of the coordinate vectors X_ u and X_ v at the 
edge of the patch are required then atleast a fourth order PDE is needed to generate the surface. The 
PDE may be written as 


\&_ , ijn 2 

du 2 dv 2 


X = 0 


(2.5) 


where X = (x(u,v), y(u,v), z(u,v)). 

The appropriate boundary conditions for Eq. (2.5) is the value of X and its normal derivative around 
the edges of the domain in the (u,v) plane. Since the generating Eq. (2.5) is an elliptic PDE, the solution 
becomes very sensitive to the choice of boundary conditions. The boundary conditions act as a powerful 
tool for surface manipulation by a designer and can be used as a design parameter in an optimization 
process. The boundary conditions on function X_ are so chosen that the curves forming the edges of the 
surface patch have the designed shape. The directions of the vector 2L U and 2L V are tangential to the 
isoparametric lines on the surface. Therefore by altering the values specified for X_ u and X_ v along the 
boundaries one can effect the direction in which the surface moves away from the edges of the patch. 
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Fig. 2.1 Six control points NURBS wing section. Fig. 2.2 Effect of control point movement. 



Fig. 2.3 Effect of increasing the weight of a control point 



Fig.2.4 Control point mesh for an M6 wing. Fig. 2.5 NURBS surface for an M6 wing. 
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The general solution of Eq. ( 2 . 5 ) can be written in the form 


OO 


X = A 0 (u) + 5^(A n (u) cos(nv) + B n (u)sin(anv)) 

n = l 

where the coefficient function A„(u) and B„(v) are of the form 


A n (u) = a nl e anu + a n2 ue anu + a n3 e- anu + a n4 ue- anu 


B.(v) = b nl e nv + b n 2 ve nv + b n 3 e-" v + b n 4 ve-“ v 

and a n]> an2, a n 3, a^, b n i, b n2 , b n 3, b„4 are vector valued constants that can be found for a particular 
solution by Fourier analysis of the condition imposed on the isoparametric lines bounding the patch. 
Now consider the generation of PDE surfaces for the HSCT configuration. It is necessary to set up the 
problem as a boundary value problem in (u,v) space, with boundary condition specified along curves in 
the (u,v) plane. One of the boundary curve is taken to be the plan outline of the intersection of the wing 
and the fuselage. On this curve, u is taken to be zero and the shape is given parametrically in terms 
of v. Another boundary curve is taken to be u=l and again is given parametrically in terms of v. This 
curve is an airfoil section at the mid portion of the double delta wing. The whole wing is generated as 
two separate sections and the fuselage is represented as a Fourier series 

Families of wing sections are described by combining a mean line and a thickness distribution. The 
resultant expression possesses the necessary features that suit the problem, mainly the concise description 
of a wing section in terms of several design parameters. The design parameters are: M = the maximum 
ordinate of the mean line or camber, and C = chordwise position of maximum ordinate. The thickness 
distribution is a bit different from the regular NACA four-digit wing section representation. Here the 
design parameters are T = The maximum thickness, PI = first Fourier wing shape parameter, and 
P 2 = second Fourier wing shape parameter. The £ - coordinate is first mapped into the chord line 
x = x(r) = x(MO) forward and the reversed to cover both the top and bottom of the section. The 
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mean line equation is 


y-c{x) = ^{2Cx-x*), x <c 

(1 - 2C + 2Cx - X s ) 
y c (x) = M-i x > C 


( 2 . 6 ) 


(i-cy 


The section thickness is given by 


y T (x) = (sin(2x) + Pi sin(4x) + P 2 sin(6x)) 


(2.7) 


The section coordinates are 

xi (r, Pi) = x , y! (f,Pl) = y c (x) ± yr(x) 

The PDE surface generated is shown in Fig. 2.6 and the various parameters that define the surface 
are shown in Fig. 2.7. 

23 Graphic Interface 

To get a feel of how the PDE surface behaves with the change in design variables, a graphic interface 
software Airplane Design and Shape Parametrization (ADSP) is developed. The software is based on 
Forms library and runs on a Silicon Graphics machine. The design variables are represented in the form 
of buttons for the convenience of user. A user can select anyone of these design variables and change 
to see interactively the surface being changed. A separate window is provided to show the values of 
the design variables being changed. The program is also menu driven and the object can be represented 
either as shaded polygon or as wire frame. Figure 2.8 shows a snapshot of the software. 

2.4 Unstructured Grid Generation and Method of Solution 

VGRID3D [10] advancing front grid generation code is used to generate unstructured grids. The 
procedure begins by inputting points defining the surface of the geometry. The surfaces are decomposed 
into smaller patches (3/4-sided) and then a background grid that defines the local grid characteristics 
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Fig. 2.6 PDE surface mesh for HSCT type configuration. 



Fig. 2.7 Surface parametrization of HSCT type configuration. 
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Fig. 2.8 Snapshot of the PDE graphic interface software 
"Airplane Design and Shape Parametrization" 


ll 



such as grid spacing, stretching and stretching direction is specified. The advancing front algorithm first 
places points on the boundary segments that define the solution domain. This yields the initial front. 
Using the information stored on the background grid, the surface patches are first triangulated and this 
forms the interior surface triangulation. The triangulation is not on the actual NURBS or PDE surfaces, 
but it is close to the surface. The resulting triangulation is then projected on to the actual surface using 
the GRIDTOOL [13] which is then used to generate the volume grid. 

The governing equations are the three-dimensional unsteady Euler equations for inviscid compressible 
flow. For a bounded domain 9 with a boundary 09, the time dependent Euler equations in integral form 
can be written as : 


where 


lTtjjj® dv + If = 0 

n an 

Q = {p,pu,pv,pw,pe 0 } T 
and 

F(Q) = 


( 2 . 8 ) 


' (pu,pv,pw) \ 

[pu s + P , puv , puw) 

< (pvu,pv 2 + P,pvw) > 

(pwu, pwv, pw 2 + P) I 

. ((pe 0 + P)u,(pe 0 + P)v,(pe 0 + P)w) ) 

In the preceding equations, p is the density, u, v and to are the x, y, z components of the velocity, e 0 
is the total energy per unit volume and P is the pressure. The equations are nondimensionalized by a 
reference density poo and speed of sound a <». Assuming an ideal gas, the pressure is written as 


P = (7 - 1) 


^e 0 - \p ( u2 + y2 + ™ 2 )) 


where 7 represents the ratio of specific heats. 


(2.9) 
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The inviscid flow field is computed on the unstructured grids using USM3D, a three-dimensional 
upwind flow solver developed at NASA/LaRC [14]. The spatial discretization is accomplished with a cell 
centered finite volume formulation using the flux difference splitting procedure. The solution is advanced 
in time using a three-stage Runge-Kutta time stepping scheme. Local time stepping and implicit residual 
smoothing are used to accelerate the convergence of the solution to a steady state. 

Unstructured grid is generated around the NURBS M6 wing using the advancing front technique 

O 

and is shown in Fig. 2.9. Converged solution is obtained for Moo=0.84 and a=3.06 and the upper 
surface pressure contour is shown in Fig. 2.10. The figure clearly shows a double shock wave on the 
upper surface [15]. 

The surface triangulation for the HSCT type PDE surface is shown in Fig. 2.11. To simulate this 
HSCT configurations with engines and to study the performance features, two nozzles with square cross- 
section are placed just below the wings. The surface triangulation for this configuration is shown in Fig. 
2.12. Both configurations are tested for Moo =0.84 <tnd <^—5 . Three stage Runge-Kutta time stepping 
scheme is used to obtain a converged solution. Figure 2.13 shows the shaded Cp plot for the HSCT 
configuration without engines. Contours are plotted by taking a cutting plane at the mid section of the 
configuration. A shock wave is seen at the upper surface of the wing. A total lift of 0.33358 and a drag 
0.04301 were obtained. Figure 2.14 shows the Cp plot for the HSCT with engines. Similar results were 
obtained, but the lift was found to be 0.313434 and the drag 0.05932. 
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Fig. 2.11 Unstructured grid for the HSCT type configuration. 



Fig. 2.12 Unstructured grid for the HSCT type configuration 
with engines. 



Fig. 2.13 Cp plot for the HSCT configuration without 
engines. 



Fig. 2.14 Shaded Cp plot for the HSCT type 
configuration with engines. 
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3. GOVERNING EQUATIONS FOR 
POTENTIAL FLOW SOLUTION 


A low-order potential-flow panel code for modeling complex three-dimensional geometries is used 
to calculate surface pressure variations. The flow field is assumed to be inviscid, irrotational and 
incompressible. The velocity potential is given by the Laplace’s equation: 

V 2 # = 0 (3.1) 

The potential at any point P may be evaluated by applying Green’s Theorem which results in the following 
integral equation 

J (* - *i)n. JJ n(y* - V*,)^ (3-2) 

S+W + Sn S+W+S„ 

It is assumed that the wake is thin and there is no entrainment, so the source term for the wake disappears 
and the jump in normal velocity across the wake is zero. Hence the simplified equation becomes 

* P = hJJ ( * " V (?) dS ~HJ (?) V$ " V * ,)d5 

s /i\ m (33) 

+ i// + 

w 

Dirichlet type boundary condition is used to solve Eq. (3.3). The total potential $ can be viewed as 
being made up of an onset potential $oo and a perturbation potential <f> = <P - #oo- The potential of 
the fictitious flow is set equal to the onset potential, 4>oo- With this boundary condition, the singularities 
on the surface tend to be smaller than if the potential of the fictitious flow is set to zero because the 
singularities only have to provide the perturbation potential instead of the total potential. The general 
equation for the potential at any point P cab be written as 

$p =[ JJ n n . y + Kn v ] + //GW/ Hwn- V (?)" "f" ^oop (3*4) 

S-P s w 


17 



where K = 0 if P is not on the surface, K = 2 tt if P is on a smooth part of the outer surface, and 
K = - 2 7r if P is on a smooth part of the inner surface. If the surface is broken up into panels, Eq. 
(3.4) can be written in discretized form, breaking the integrals up into surface integrals over each panel. 
A constant strength source and doublet distribution is assumed over each panel and so the doublet and 
source strengths are factored out of the integrals. Taking point P to be at the centroid on the inside of one 
of the panels, the surface integrals over each panel are summed for all panels. For the panel containing 
point P, the surface integral is zero and only the -2x[i T term remains in the bracketed part of Eq. (3.4). 
For all other panels, the surface integral is used and the — term is zero since the point P is not 
on the surface of any other panels. The process is repeated for point P at the centroid of every panel to 
yield a set of linear simultaneous equations to be solved for the unknown doublet strength on each panel. 
The surface integrals represent the velocity potential influence coefficients per unit singularity strength 
for panel K acting on the control point of panel J. Hence Eq. (3.4) becomes 

N, N B Nw 

T: {hkCjk) + ^ {vkBjk) + ^2 ((J-w l Cjl) = 0 (3-5) 

K = 1 K = 1 L= 1 


where 

Bjk = jj\dS (3.6) 

K 

and 

K v ’ 

Cjj = -27T 

The coefficients Cjk and Bjk represent the velocity potential influence coefficients per unit singularity 
strength for panel K acting on the control point of panel J. Equations (3.6) and (3.7) are functions of 
geometry only and thus can be solved for all panels to form the influence coefficient matrix. Since the 
source values are known, they may be transferred to the right hand side of the matrix equation. Solutions 
for Eqs. (3.6) and (3.7) can be found in [15], 
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4. COMPUTATION OF SENSITIVITY TERMS 
WITH RESPECT TO DESIGN PARAMETERS 


After discretization, the steady state governing equations of the fluid flow and the boundary conditions 
can be expressed as 

R(Q(D),X(D), D) = 0 (4.1) 

where Q is the velocity potential for a potential flow and X is the computational grid or panels in the 
present case, and D is a vector of independent input (design) variables. 

As an alternative approach a quasi-Newton iteration can be applied which can be represented as 

-|^AQ = R“ (4.2) 

Q» + l _ QU + A Q (4.3) 

Application of the quasianalytical methods that has been investigated by many researchers [16] 
requires the construction and evaluation of many derivatives (eg. the Jacobian matrices and §y). 
For advanced CFD codes the task of constructing exactly all of these required derivatives by hand and 
then building the software for evaluating these terms is extremely complex, error prone, and practically 
speaking, impossible. A promising possible solution to this problem may be found in the use of a 
technique known as automatic differentiation (AD) which has been used in the present study. The process 
involves the application of a precompiler software tool that automatically differentiates the application 
program source code from which sensitivity derivatives are to be obtained. The output of the AD 
precompiler procedure is a new source code which, upon compilation and execution, will compute the 
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numerical values of the derivatives of any specified output functions with respect to any specified input 
parameters. This AD precompiler tool has been efficiently tested by Bischof et al. [17] and Green et al. 
[18] in applications to an advanced CFD flow-analysis code called TLNS3D [19]. 

When AD is applied directly to the potential flow code, the resulting AD-enhanced code calculates 
the required sensitivity derivatives through an iterative process. From the discussion in [17], the process 
whereby sensitivity derivatives are iteratively calculated after the application of AD can be represented 
conceptually by combining Eq. (4.2) and (4.3) i.e., the basic CFD solution procedure, and differentiating 
with respect to D. The result is 


Q m+1 = Q ,n - P"R' n - P ,n R“ 


where P = 


( 3R 

ISQ 


)-• 


The Automatic Differentiator (ADIFOR) procedure generates a new version of the potential flow 
code that has the capability to calculate the derivatives of lift, drag, and pitching moment with respect 
to a wide variety of different types of input parameters (including parameters related to the geometric 
design). Table 4.1 show the non-geometric and geometric ADIFOR sensitivity values compared with 
finite difference. It is seen that the results obtained by ADIFOR is in good agreement with that of the 
finite difference. 


4.1 Grid Sensitivity with respect to Design Parameters 

Typical CFD calculations are performed on a computational mesh that is “body-oriented”. Changes 
in the geometric shape result in the movement of grid points throughout the entire mesh. One method for 
calculating these grid sensitivity terms is by finite differences. If forward difference approximations are 
selected, for example, the mesh generation code is used to produce one additional perturbed grid for a 
slightly perturbed value of each geometric shape design variable of interest. This procedure is generally 
expected to be reliable in producing accurate grid sensitivity terms because the relationships that are 
associated with the mesh generation process should be very smooth by design. 
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In this study, however, the precompiler tool AD (i.e., ADIFOR) is applied directly to the grid 
generation program to successfully calculate the grid sensitivity terms. Figures 4.1a and 4.1b show the 
comparison of the grid sensitivity with respect to camber between the ADIFOR and finite difference 
results. These grid derivatives were subsequently coupled directly to the AD-enhanced potential flow 
code PMARC. The final result is the successful calculation of aerodynamic sensitivity derivatives with 
respect to geometric design parameters. 

4.2 Optimization Problem 

An objective of multidisciplinary optimization of a vehicle design is to extremize a payoff function 
combining dependent parameters from several disciplines. Most optimization techniques require the 
sensitivity of the payoff function with respect to free parameters of the system. For a fixed grid and 
solution conditions, the only free parameter are the surface design parameters. Therefore, the sensitivity 
of the payoff function with respect to design parameters are needed. 

The present optimization strategy is based on maximizing the lift coefficient. Cl, in response to 
surface perturbation, subject to pre-determined design constraints. Upper and lower bounds are set for 
each design parameter and the sensitivity derivatives of the objective function, and the constraint, 
are obtained as previously described. Throughout the analysis, the drag coefficient, Cd. is to be 
no greater than the value of the initial design. The strategy, illustrated in Fig. 4.2, requires that the grid 
and grid sensitivity derivatives be provided dynamically during the automated optimization process. 

Optimization on the HSCT type configuration shown in Fig 2.6 was carried out on a SGI machine 
with a memory capacity of 512 MB. Sixteen design variables were selected for the optimization process. 
A total of twelve design optimization cycle was performed and each iteration took approximately 7.5 
min of cpu time. It was seen that the lift which was initially 0.0356 became 0.1245. The initial and final 
shapes with shaded Cp plot are shown in Figs. 4.3a - 4.3c. 
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5. CONCLUSIONS 


A feasibility study has been conducted for using non-uniform rational B-splines and partial differential 
equation to represent aerospace vehicle. Unstructured grids generated by advancing front method is used 
to obtain surface triangulation and converged solution obtained using an Euler solver. It is seen that 
the surfaces produced by PDE tend to be very smooth and could accurately represent the blending of 
one component into another (wing into fuselage). The NURBS surface had the advantage of altering 
the surface by the movement of the control points. It is also noted that the PDE method could easily 
generate a wide variety of aerospace configurations with a slight change in the design parameter and 
hence is very useful in an optimization process. 
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Table 1. Comparison Of ADIFOR results with finite difference for geometric 
and nongeometric design variables. 



Fig. 4.1a ADIFOR Y-coordinate Fig. 4.1b Finite diff. Y-coordinate 

sensitivity with respect to camber. sensitivity with respect to camber. 
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Fig. 4.3a Potential flow on the original pig, 4 . 3 b Potential flow on the optimized 

configuration. configuration. 



Fig. 4.3c Shaded plot of the original and optimized HSCT type 
Configuration. 
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